arXiv: 1506.0632lv2 [quant-ph] 18 Nov 2015 


Qubit measurement error from coupling with a detuned neighbor in circuit QED 


Mostafa Khezri, 1,2 ’* Justin Dressel, 1 and Alexander N. Korotkov 1 

1 Department of Electrical and Computer Engineering, 

University of California, Riverside, CA 92521, USA. 

2 Department of Physics, University of California, Riverside, CA 92521, USA. 

(Dated: November 19, 2015) 

In modern circuit QED architectures, superconducting transmon qubits are measured via the 
state-dependent phase and amplitude shift of a microwave field leaking from a coupled resonator. 
Determining this shift requires integrating the field quadratures for a nonzero duration, which can 
permit unwanted concurrent evolution. Here we investigate such dynamical degradation of the 
measurement fidelity caused by a detuned neighboring qubit. We find that in realistic parameter 
regimes, where the qubit ensenrble-dephasing rate is slower than the qubit-qubit detuning, the joint 
qubit-qubit eigenstates are better discriminated by measurement than the bare states. Furthermore, 
we show that when the resonator leaks much more slowly than the qubit-qubit detuning, the mea¬ 
surement tracks the joint eigenstates nearly adiabatically. However, the measurement process also 
causes rare quantum jumps between the eigenstates. The rate of these jumps becomes significant 
if the resonator decay is comparable to or faster than the qubit-qubit detuning, thus significantly 
degrading the measurement fidelity in a manner reminiscent of energy relaxation processes. 

PACS numbers: 03.67.Lx, 85.25.-j 


I. INTRODUCTION 

Recent years have witnessed the rapid evolution of su¬ 
perconducting circuit QED technology for quantum com¬ 
putation [1-9] (reviewed in [10-12]). The most recent de¬ 
velopments have converged on charge-insensitive designs, 
based on transmons [13], which can be dispersively mea¬ 
sured with coupled microwave resonators. Multi-qubit 
chips based on these designs have recently demonstrated 
high-fidelity entangling gates [4, 5], which are now near¬ 
ing the gate fidelity thresholds necessary for implement¬ 
ing practical quantum error correction protocols [14-16]. 
Indeed, several groups have recently demonstrated bit- 
flip error correction in such multi-qubit superconduct¬ 
ing processors [17-20]. With the gate fidelity reaching 
such unprecedented levels, it is now interesting to iden¬ 
tify and address more subtle sources of error that can 
arise in such a multi-qubit environment, such as the ef¬ 
fect of non-tunable qubit-qubit or qubit-bus coupling [21] 
on the dispersive measurement fidelity. 

Unlike the textbook projective measurements usually 
assumed in the quantum computing literature [14] , which 
involve instantaneous state collapse, realistic measure¬ 
ments occur over a nonzero duration of time. In the 
transmon-based circuits we consider here, each qubit 
is dispersively coupled to a pumped microwave res¬ 
onator such that the leaked field is phase-shifted (and, in 
general, amplitude-shifted) by a qubit-state-dependent 
amount [1]. The leaked field is then passed through 
an amplifier and mixed with a local oscillator to pro¬ 
duce a noisy homodyne signal. This signal needs to be 
integrated until the signal-to-noise ratio exceeds an ac¬ 
ceptable discrimination threshold (more advanced signal 
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processing techniques can moderately increase the mea- 
surement fidelity [22]). For an isolated qubit, the in¬ 
crease of the signal-to-noise ratio by longer integration 
is limited by the energy relaxation (and excitation) pro¬ 
cesses. However, in circuits intended for quantum com¬ 
putation, the qubits will also be coupled to frequency- 
detuned neighbors, which may permit unwanted dynam¬ 
ics to additionally degrade the measurement fidelity. We 
wish to better understand the detailed dynamics of a re¬ 
alistic transmon qubit measurement, and minimize the 
dynamical measurement error that will arise from the 
coupling to a neighboring qubit (or bus, which plays a 
similar role). 

In this paper, we demonstrate that in typical exper¬ 
imental parameter regimes, where the qubit ensenrble- 
dephasing rate due to measurement is slower than the 
qubit-qubit detuning, dynamical measurement error al¬ 
ways exists when distinguishing the bare energy states of 
coupled qubits. However, this measurement error can be 
decreased by distinguishing not the bare energy states, 
but instead the qubit-qubit eigenstates that are station¬ 
ary under the effect of the qubit-qubit coupling and de¬ 
tuning (similarly to the measurement of coupled phase 
qubits analyzed in [21]). Despite the fact that only the 
main qubit is being measured, the relatively slow mea¬ 
surement process allows the two-qubit system to collapse 
to these stationary eigenstates, in contrast to what may 
be naively expected from textbook projective measure¬ 
ments. Notably, these eigenstates have also been shown 
to be a natural choice for the logical encoding of high- 
fidelity multi-qubit gates [21] (for similar reasons), which 
makes multi-qubit eigenstates an unambiguously optimal 
choice for logical encoding in realistic parameter regimes. 

We further demonstrate for coupled transmon mea¬ 
surement that, in addition to the ensenrble-dephasing 
rate and the qubit-qubit detuning, the measurement fi- 



2 


delity depends on a third important parameter: the read¬ 
out resonator energy decay rate due to leakage into a 
transmission line. For decay rates much slower than the 
qubit-qubit detuning (as is typical in experiments, e.g., 
[18-20]), the leaked resonator field nearly adiabatically 
follows the qubit-qubit eigenstate to produce little er¬ 
ror. However, for decay rates that are comparable to or 
larger than the qubit-qubit detuning, the resonator de¬ 
cays more rapidly than it can equilibrate with the qubit- 
qubit eigenstates, causing frustrated dynamics during the 
measurement. Such rapid resonator decay will primarily 
couple the leaked field to the bare energy states, while 
the fast inter-qubit oscillations (compared to the mea¬ 
surement rate) will relate the output signal to the joint 
qubit-qubit eigenstates. This frustrated dynamics leads 
to random quantum jumps between the eigenstates. We 
derive the rate of these quantum jumps (which we call 
a switching rate) using a model based on fluctuations of 
the photon number in the resonator, which perturb the 
two-qubit eigenstates. We show that the switching can be 
significant for rapid resonator decay, but becomes almost 
negligible for realistically slow resonator decay. We also 
derive the measurement error probability resulting from 
these quantum jumps, and show that it accumulates al¬ 
most linearly with integration time in an analogous way 
to the error from energy-decay (T)) processes. 

This paper is organized as follows. In Section II, we 
introduce the considered system, formulate the problem, 
and discuss how to model the ensemble-averaged dynam¬ 
ics. In Section III, we identify three qualitatively distinct 
parameter regimes in the ensemble-averaged dynamics: 
textbook, adiabatic, and frustrated. In Section IV, we 
study the transition between the adiabatic and frustrated 
regimes as the resonator decay is varied, by introducing 
a simple model of a semiclassically fluctuating field in 
the resonator that produces random quantum jumps be¬ 
tween the eigenstates. We derive the average switching 
rate for these jumps, and numerically confirm this jump 
behavior by simulating quantum trajectories in the fast 
resonator decay regime. In Section V, we demonstrate 
that the contribution of these jumps to the measurement 
error is nearly linearly increasing with integration time, 
and find the error minimized over the integration time. 
We conclude in Section VI. 


II. CONSIDERED SYSTEM AND ITS 
ENSEMBLE-AVERAGED EVOLUTION 

The archetypal circuit QED system we consider here 
is shown in Fig. 1. A superconducting transmon (main 
qubit) with frequency w q (which includes the ac Stark 
shift) is capacitively coupled to a driven readout res¬ 
onator with bare frequency w r , and is also capacitively 
coupled to another transmon (neighboring qubit) with a 
detuned bare frequency w n , such that the qubit-qubit de¬ 
tuning A = w q — w n is much larger than the qubit-qubit 
coupling g : |A| g. (The role of the neighboring qubit 
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FIG. 1. (Color online) Analyzed system. A measured (main) 
transmon qubit (blue) with frequency ta q has capacitive cou¬ 
pling jtoa detuned neighboring qubit (red) with a frequency 
a) n such that g -C |A|, where A = to q — ui n (io q includes the ac 
Stark shift). A readout resonator with frequency co r is disper- 
sively coupled to the main qubit, and thus is frequency-shifted 
by ±x depending on the main qubit state. During measure¬ 
ment, this resonator is driven with a coherent microwave e 
at a frequency a>d- The field leaks from the resonator (with 
the energy decay rate k) to a transmission line, where it is 
amplified and mixed with a local oscillator to measure the 
quadrature I ( t ) that is sensitive to the qubit-state-dependent 
phase and amplitude shift. The coupling of the main qubit 
with the neighboring qubit contributes to the qubit measure¬ 
ment error. 


can be played by a bus resonator; we consider a qubit for 
definiteness.) We assume that the Purcell decay [23-25] 
of the main qubit through the resonator has been sup¬ 
pressed by a filter [26-28], and that the resonator and 
main qubit are sufficiently detuned to treat their coupling 
as effectively dispersive (implying the rotating wave ap¬ 
proximation) [1, 29]. We also assume that the transmon 
energy levels outside of the qubit subspace are taken into 
account through renormalization of the state-dependent 
dispersive shift ±y of the resonator frequency. The read¬ 
out resonator is additionally driven by a coherent field e 
at a microwave frequency Wd, which then leaks to a trans¬ 
mission line at an energy-decay rate k (the setup can be 
either in transmission or in reflection). The leaked field 
is passed through an amplifier and mixed with a local 
oscillator to perform a homodyne measurement, which 
isolates the qubit-state-dependent phase and amplitude 
shift caused by the dispersive coupling; the information¬ 
carrying quadrature is denoted as I(t) in Fig. 1. 

While this measurement procedure is largely under¬ 
stood for a single qubit coupled to the readout resonator 
[1, 29-31], we investigate here how the addition of the 
neighboring qubit will contribute to the measurement 
error. Specifically, we wish to find out whether the 
wavefunction “tail” probability {g/ A) 2 contributes to the 
measurement error or not. In this paper we focus on 
discriminating the bare states 110) and |00) (with qubit 
ordering convention [main, neighbor)) or the states 110) 
and |00), where 110) is the eigenstate that accounts for 
the qubit-qubit interaction. We assume the logic state 
of zero for the neighboring qubit for simplicity, without 
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significant loss of generality, because the discrimination 
of the states 111) and |01) (or 111) and |01)) is a very sim¬ 
ilar problem. Also, the discrimination of all four states 
in the case when both qubits are measured is a simple 
generalization of our basic problem. Note that we do not 
consider another important question: deterioration of a 
superposition a|00) + /3|01) (or a|00) + /?|01)) after the 
main qubit measurement; however, the mechanism of this 
deterioration is similar to what we consider. Also note 
that in an architecture [21, 32], in which the zero state 
of the neighboring qubit is used as a resource to decrease 
crosstalk, our assumption of discriminating 110) and 100) 
(or 110) and |00)) is naturally satisfied. We will refer to 
the pair of states to be discriminated as logical 1 and 0. 

The logical states are discriminated by integrating the 
fluctuating output signal I ( t ) over time and then compar¬ 
ing the result with a threshold. Therefore, the discrimi¬ 
nation error P err (discussed in more detail in Sec. V) can 
be calculated from the “overlap” of the probability distri¬ 
butions of the integrated result for the two logical states. 
The error depends on the chosen threshold (in Sec. V we 
will consider the symmetric and optimal thresholds) and 
on the integration time. As will be discussed later, the 
measurement error P err has a minimum as a function of 
the integration time, which is determined by the rate of 
“switching” (quantum jumps) between the qubit states, 
resembling the energy relaxation events. 

Our final goal is to find such optimized measurement 
error for distinguishing the bare-basis states 110) and 
100), and for distinguishing the eigenbasis states 110) 
and 100), thus finding which encoding basis is prefer¬ 
able in the circuit QED measurement. The analysis of 
a similar question for the measurement of phase qubits 
showed [21] that using the eigenbasis is preferable. In 
this paper we will obtain a similar result for the parame¬ 
ter regime of typical circuit QED measurements [18-20], 
even though the measurement dynamics is significantly 
more complicated than for phase qubits. In particular, 
we will show that in contrast to what is expected for 
a textbook projective measurement, the bare-basis er¬ 
ror exceeds (g/ A) 2 /2, while there is no such limitation 
for discriminating the eigenstates. For the eigenbasis the 
limitation comes from the quantum jumps between the 
eigenstates 110) and |01); however, for typical experimen¬ 
tal parameters this limitation is almost negligible. 

To obtain these results, we first discuss how to model 
both the coherent and incoherent aspects of the evolution 
for the ensemble-averaged case. This ensemble-averaged 
dynamics will be sufficient to identify broad parameter 
regimes of interest for the coupled-qubit measurement, 
and to identify which qubit-qubit bases are preserved by 
the measurement in these regimes, but will be insufficient 
for understanding and quantifying the measurement error 
for specific realizations. In Section IV, we will generalize 
the ensemble-averaged approach to model the individ¬ 
ual quantum trajectories, which will allow us to under¬ 
stand and derive the measurement error induced by the 
qubit-qubit coupling. Note that we consider only one 


neighboring qubit, while in practical architectures (e.g., 
in surface codes) there are several neighbors; however, 
the generalization of our theory to several neighboring 
qubits is rather straightforward. 


A. Coherent evolution 


The total resonator-qubit-qubit Hamiltonian can be 
split into five terms, 


H = H r +H q + H d + H ql + H qq . 
The bare-energy contributions (H = 1) are 


H r = w r a^a, 


H _ ZV!) + A1V2) 
— 2 uz ' 2 ° z ’ 


(1) 

( 2 ) 


where aP = 11)(1|y — |0)(0|j are the Pauli z operators 
for each qubit (j = 1,2), a* (a) are the raising (lower¬ 
ing) operators of the coupled resonator mode that satisfy 
[a, a'] = 1, and cu q , u ni and w r are the bare frequencies of 
the main qubit, neighboring qubit, and resonator. The 
resonator drive contribution has the form 


H d (i) = e e -i “ dt a* + £* e iuidt a. (3) 


The bare energies are modified by the dispersive qubit- 
resonator coupling 

H v = xvP^a, (4) 

which shifts the resonator frequency by ±y depending 
on the qubit state or, alternatively, shifts the qubit fre¬ 
quency, 

W q = UJ q + SCO q , (5) 

by the ac Stark shift 6ui q , depending on the number of 
photons in the resonator (we include the Lamb shift [1, 
33, 34] into u q ). The qubit-qubit coupling Hamiltonian 
(assuming the rotating wave approximation) is 

F qq = 5 (|01)(10| + |10)(01|), (6) 

and we are interested in the case of strongly detuned 
qubits, g <C |w q — w n | (for simplicity we assume g > 0). 

Note that the qubit-qubit coupling in Eq. (6) coher¬ 
ently mixes the single-excitation subspace {101), 110)} 
and produces the eigenstates of the qubit-qubit Hamilto¬ 
nian H q + H qq that are rotated from the bare states by 
an angle 0=1 arctan(2g/A), 


110) = cos0|lO) + sin0|Ol) 
|01) = cos0|Ol) — sin0|lO) 
where 



A = W q - W n 


( 8 ) 
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is the (ac Stark-shifted) qubit-qubit detuning and the 
approximation is to lowest order in <//|A| <C 1. If the 
measurement process occurs effectively in this eigenbasis, 
then an initially bare state 110) will collapse into the in¬ 
correct eigenstate 101) with an error probability (g/ A) 2 , 
resulting in additional measurement error. In Section III 
we will clarify which parameter regimes of the measure¬ 
ment naturally select the eigenstates of Eq. (7) in this 
manner. Note that for brevity of notations, in inequali¬ 
ties describing the parameter regimes we will use A in¬ 
stead of |A|. 


B. Incoherent evolution 


that describes the ensemble-averaged evolution for a du¬ 
ration d t. The Hamiltonian evolution is then interleaved 
between these incoherent updates: p' >->• Up'W. Ex¬ 
panding the full increment dp = Up'U' — p to linear 
order in dt produces the standard master equation form 
of Eq. (9). We note, however, that the update in Eq. (10) 
is not only conceptually transparent, but (as we checked) 
is more numerically stable for simulation purposes and in 
some regimes is faster than solving Eq. (9). 

If additional decay channels are present, they can 
be added phenomenologically to the incoherent sum in 
Eq. (10). For example, qubit energy-decay with rate l/T) 
and environmental qubit dephasing with rate T e have the 
forms 


In addition to the coherent evolution given by the 
Hamiltonian in Eq. (1), the energy in the resonator in¬ 
coherently decays to a transmission line at the rate k. 
Assuming that all leaked photons may not later return 
to the resonator, we can model the ensemble-averaged 
Markovian evolution of the joint qubit-qubit-resonator 
state with a master equation [33] 


dp = — i[H , p]dt + (ndt) apa) 


ndt 

~Y~ 


{a'ap + pa^a). (9) 


Physically, we can interpret this equation as stating that 
in a small interval df the system does not only coherently 
evolve with the usual evolution operator 


U = exp(— iHdt), 


but additionally has one of two distinct incoherent pro¬ 
cesses happen (e.g., Ref. [35]): 

(a) Each of N photons in the resonator may escape with 
probability ndt 7 which modifies the resonator state 
with the decay operator 

Afdecay = VK~df d. 


M T1 = \jdt/T[ (a) (b) * * * * * * * * * 1] <T ( y + yjdt/T™ a {2 \ 

Mdephase = \] df O ^ + \/l1, 2) df (J (2) , 

which will modify the null result operator accordingly to 
include all decay channels Mk 

M nu „ = y/l - EMM,, 

These additional decay channels correspondingly modify 
the linear increment in Eq. (9) in the standard way. For 
simplicity, we will neglect such additional decay channels 
in most of what follows, in order to focus solely upon 
the effects of the neighboring qubit on the measurement 
fidelity. When we do add these effects, we will assume 
that t[ 1] = T[ 2) = Tj and T^ = T^ 2) = T e . 

Now let us briefly review some results for measure¬ 
ment of a single qubit [1, 29-31], which we will use as 
a starting point and to introduce notations. For the 
qubit in the state |1) or |0), the effective frequency of 
the resonator is w r ± \ (the upper sign is for the state 
11)). Then the evolution of the resonator coherent state 
|a±) = e _ l“ ± l 2 / 2 o 1 ^.(n\)~ 1 / 2 e~ muidt \n) (we use the 

rotating frame e~ lWdt ) is 


(b) All N photons stay in the resonator with probability 
1 — ftdfiV, which modifies the resonator state with 
the null result (no decay) operator 

M nu n = i/l — ft df cda. 


d± = -i{ A r ± x) a±-^a±- is, (H) 

where A r = w r — Wd is the bare resonator-drive detuning. 
The steady-state solution of this equation is 


These measurement (Kraus) operators for the incoher¬ 

ent part of the evolution satisfy the usual completeness 

condition M] ecay M decay -(-M^ ull M null = 1, indicating that 

the probabilities for each possibility to occur are correctly 

normalized [14]. Mixing together both possibilities (i.e., 

discarding any record of whether the decay happened or 

not) produces the updated mixed state 

p' = Miecay p M \ ecay + M null pM^ uU , (10) 

= k d t apa) + x/u — K d tad a p x/l — k d tad a. 


a± = 


— l£ 


k /2 + *(A r ± x) ’ 

and the corresponding mean photon number is 


( 12 ) 


n± = a±| =n n 


k 2 + 4(A r ± x ) 2 ’ 


(13) 


which we expressed via the photon number at exact res¬ 
onance, h max = 4|£| 2 /k 2 . The ac Stark shift is then [29] 


Sui q = 2y Re(a^_o:_) 


(14) 
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and the measurement-induced ensemble dephasing rate 
is [29] 


r m = 2xlm(a+a_) = re ^ ^ . (15) 

These results can be expressed in terms of n± and n max 
as 


5cu q = 2% — 




4(A? - X 2 ) 


r m = 


^max 

8x 2 n+n- 


(16) 

(17) 


which reduce to the simple formulas [1, 29] dwq ss 2\n, 
r m ~ 8 xrn/n when n + « ? 1 _ ss h max . One of the ways 
to interpret the measurement-induced dephasing process 
is as being caused by fluctuations of the ac Stark shift 
that arise from the fluctuating photon number. The to¬ 
tal ensemble-depliasing rate T = T m + T e generally in¬ 
cludes additional environmental dephasing T e , but we 
will mostly neglect T e for simplicity. The measurement- 
induced ensemble dephasing is related to the distin- 
guishability time (sometimes called the “measurement 
time”) 


r=(27 7 r m )- 1 , (18) 

needed for achieving unit signal-to-noise ratio in the 
quadrature output, where 77 £ [ 0 , 1 ] is the quantum ef¬ 
ficiency of the detection circuit. 

We emphasise that these standard results for <5w q and 
T m are for the measurement of a single qubit; moreover, 
they implicitly assume the “bad cavity limit” in the sense 
that the qubit evolution is much slower than re (in this 
case it is sufficient to consider only coherent states in 
the resonator, entangled with the qubit, which leads to 
relatively simple formulas). Therefore, we should not ex¬ 
pect that these results are directly applicable to our prob¬ 
lem, which focuses on evolution involving the neighboring 
qubit. In particular, when the qubit-qubit detuning A is 
larger than re, the relatively slow fluctuations of the pho¬ 
ton number in the resonator will not produce the same 
dephasing T m between the states 110 ) and | 01 ) (as would 
be expected for infinitely fast fluctuations). Similarly, 
for A > ft the ac Stark shift contribution to A is sup¬ 
posed to be governed mainly by n + (or h_ if the main 
state is |01)) rather than given by Eq. (16). Even though 
T m in Eq. (17) does not in general describe the ensemble 
dephasing between 110 ) and 101 ), in this paper we will 
extensively use T m defined in Eq. (17) as a notation. 


III. EIGENSTATES VS. BARE STATES 

A master equation is incapable of describing the fi¬ 
delity of the qubit measurement, even in principle, so we 
will be forced to consider the individual quantum trajec¬ 
tories in Section IV. Nevertheless, even without a more 


detailed trajectory description we can already answer the 
most basic question about the qubit measurement: does 
the ensemble-averaged evolution faithfully preserve a log¬ 
ical qubit basis? 

To answer this question, we simulate the full master 
equation in Eq. (9) [equivalently, Eq. (10) can be iter¬ 
ated] starting in either a bare state 1 10 ), or an eigen¬ 
state 1 10 ), with the resonator in an initial ground state 
for simplicity (the simulation starting in the state | 00 ) is 
trivial). When starting in 110), we calculate the evolu¬ 
tion of the bare state population P 10 , and when starting 
in 1 10 ), we calculate the eigenstate population (see 
the left and right panels in Fig. 2). If one of these pop¬ 
ulations remains very close to 1 , then we infer that the 
corresponding basis is faithfully preserved by the mea¬ 
surement dynamics. 

As shown in Fig. 2, from these simulations we identify 
three parameter regimes that have qualitatively different 
behaviors (using g <C A and T e = 0): 

(a) A (T m , re) : textbook — almost stable bare state , 

(b) (T m , k)cA: adiabatic — almost stable eigenstate , 

(c) T m « A « re : frustrated — unstable eigenstate. 

The parameters used for each of these regimes are de¬ 
tailed in the caption for Fig. 2. As expected from the 
similar analysis for measurement of phase qubits [ 21 ], for 
T m < A the measurement effectively occurs in the eigen- 
basis, while the traditional (textbook) bare-basis mea¬ 
surement requires T m A. However, transnron qubits 
have an additional important parameter that has no ana¬ 
logue in phase qubits: the resonator energy-decay rate re. 
As we will see, the relative magnitudes of re and A deter¬ 
mine the “stability” of the eigenbasis. 

In the regime (a), the resonator empties and the system 
dephases much faster than the qubit-qubit evolution, so 
the bare states 1 10 ) and | 00 ) are preserved as the optimal 
logical basis, just as we would expect from a textbook 
projective measurement. That is, our numerical simula¬ 
tion in Fig. 2(a) shows that the bare state population Pi 0 
is preserved practically at 1 , while the eigenstate popu¬ 
lation P 10 (when starting with 110)) drops by roughly 
2(g/A ) 2 during the transient (collapse) evolution. [Here 
one factor of (< 7 /A ) 2 stems from the physical collapse of 
the eigenstate to an incoherent mixture of the single¬ 
excitation bare states 1 10 ) and 1 01 ), while the second 
factor (< 7 /A ) 2 comes from plotting the eigenstate pop¬ 
ulation.] At a much longer time scale the bare-basis pop¬ 
ulation gradually decreases because non-zero g makes the 
measurement not fully projective, leading to rare transi¬ 
tions (jumps) between the states 110) and |01). Note that 
for numerical simplicity in the regime (a) we simulated 
the evolution assuming re T m , so that the qubits and 
resonator remain effectively disentangled (qubit entan¬ 
glement with the emitted field is not important for the 
master equation approach). With this approximation, we 
can simplify Eq. (9) by reducing it to a two-qubit Hilbert 
space and taking into account the interaction with the 
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FIG. 2. (Color online) Blue lines: ensemble-averaged evolu¬ 
tion of the population Pio of the bare-basis state 110) when 
the evolution starts in this state (left panels) and the popu¬ 
lation Pio of the eigenstate 110) when starting in this state 
(right panels). The dashed red lines show the initial value 1 
of all blue lines for reference. The time is normalized by the 
ensemble-dephasing rate F m due to measurement; we assume 
fixed qubit-qubit coupling and detuning with gj A = 1/10 for 
all regimes, (a) Textbook regime with A T m k, using 
directly applied qubit-dephasing of T m /A = 20 for simplicity 
(i.e., assuming k — > oo). The bare state 110) is best preserved 
by the evolution, but slowly decays at the rate 2<? 2 /T m , while 
the eigenstate population Pio additionally drops by approx¬ 
imately 2(g/A) 2 . (b) Adiabatic (experimental) regime with 
(r m ,K) <C A, using k/A = 10 _1 and T m /A = 10 -2 , set by 
assuming a weak response x/ K = 3.5 x 10 -2 and a resonator 
drive wa = oj r with power tuned to produce the steady-state 
photon number h = 10. The eigenstate j 10) is best pre¬ 
served by the evolution, but slowly decays (analogously to 
the textbook regime for Pio), while the bare population Pio 
additionally drops by 2(g/A) 2 . (c) Frustrated regime with 
r m <C A <C k, using k/A = 10 and T m /A = 10 —4 , keep¬ 
ing the same y/A and n as in the adiabatic regime. The 
bare state population Pio drops by 2(g/A) 2 compared to the 
eigenstate population Pio, and both populations show rapid 
decay. The decay rate seen in regimes (b) and (c) matches 
the analytical results for averaged incoherent quantum jumps 
between the eigenstates (see Fig. 4), an example of which is 
shown here in the bare (c) plot as the overlaid dashed yellow 
curve. 


resonator by applying the dephasing with the rate r m 
to the measured qubit. In this case the transient evo¬ 
lution occurs on the time scale F ” 1 and the population 
Pio decays with the rate 2g 2 /T m . The textbook regime 
(a) is most easy to understand and analyze. However, 
we emphasize that this regime is not realized in realistic 


experiments with transmons, in which typically A T m . 

In the adiabatic regime (b), which more closely de¬ 
scribes recent experiments [18-20], the resonator empties 
and the system dephases more slowly than the qubit- 
qubit evolution, so the eigenstates 1 10 ) and 1 00 ) are pre¬ 
served as the optimal logical basis. That is, our numer¬ 
ical simulation in Fig. 2(b) shows that the eigenstate 
population Pio is preserved at almost 1 , in contrast to 
the textbook regime, while the bare state population Pi o 
(when initially 1 ) drops by roughly 2{g/A) 2 within the 
collapse timescale (this timescale is F " 1 if k T m , while 
for k < T m everything is determined by transients). 
Again, in this drop one factor of (g/ A ) 2 comes from the 
collapse into an incoherent mixture of single-excitation 
eigenstates 110 ) and | 01 ), while the second factor comes 
from plotting the bare state population. At longer time 
scales, we also observe in Fig. 2(b) that the eigenstate 
population ^10 decays exponentially at a very slow rate. 
This occurs because of rare transitions (jumps) between 
the eigenstates 110 ) and | 01 ), discussed in more detail 
later. 

The frustrated regime (c) differs from the adiabatic 
regime (b) only by the relative magnitude of the res¬ 
onator decay k and the qubit-qubit detuning A, k A. 
Nevertheless, this regime dramatically amplifies the ex¬ 
ponential decay process observed at long times in the 
adiabatic regime (b). The rapid decay seen in Fig. 2(c) 
occurs for both bare and eigenstates, so that neither of 
these bases is good for preserving a logical state. This 
occurs because fast oscillations A (compared to dephas¬ 
ing T m ) favor the eigenbasis, while even faster decay k 
makes the outgoing photons sensitive to the bare basis. 

In both regimes (b) and (c), the system collapses to the 
eigenstates, after which the state may jump between the 
eigenstates. This behavior is evidenced in Fig. 3, showing 
the ensemble-averaged evolution in the Bloch sphere rep¬ 
resentation of the qubit-qubit single-excitation subspace. 
The ratio k/A = 1 is chosen in between the regimes 

(b) and (c). The initially bare state rapidly oscillates 
around the eigenstate axis as it spirals into this axis on 
average, indicating that the initially bare state collapses 
to an incoherent mixture of the eigenstates. After that 
the exponential decay occurs along the eigenstate axis of 
the Bloch sphere, indicating that it arises solely from a 
classical mixing process that scrambles those eigenstates. 

The physical origin of the exponential decay seen in 
regimes (b) and (c) is not apparent from examining 
the ensemble-averaged behavior of the master equation 
alone, but we shall see that this decay can be interpreted 
as arising from averaging random quantum jumps be¬ 
tween the eigenstates that occur during the continuous 
measurement process. For the remainder of this paper 
we will mostly focus on understanding the transition be¬ 
tween the adiabatic regime (b) and the frustrated regime 

(c) as k is varied with respect to A. 

We also briefly note that in principle there is a fourth 
parameter regime: (d)K«A<C F m . We do not consider 
this regime here, since in this case it is difficult to clearly 
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FIG. 3. (Color online) Qubit-qubit single-excitation 
ensemble-averaged evolution, using bare Bloch sphere coor¬ 
dinates defined by x = 110)(011 + |01)(10|, y = — i(|10}(01| — 
101)(10|) and z = |10)(10| — |01)(011, and parameters g/A = 
1/10, k/A = 1, and r m /A = 1/100 (with n = 2 and cor¬ 
respondingly x/ K — 1/40). An initially bare state 110) = 
(0, 0,1) oscillates rapidly around the tilted axis correspond¬ 
ing to the eigenbasis {110), 101)}, and approaches this axis, 
indicating the gradual collapse, which produces an incoher¬ 
ent mixture of the eigenstates. At longer times, the ensemble 
averaged state continues moving along the eigenstate axis at 
a slow rate, indicating an additional classical mixing process, 
(top) 3D plot of (x , y, z) evolution, showing the spiraling evo¬ 
lution to the eigenstate axis and then along the axis, simulated 
for r m t € [0,40]. (bottom) Slice of x-z plane, with sphere sur¬ 
face shown as the dashed gray curve and the eigenstate axis 
shown as the dashed red line tilted from the bare z axis by 
the angle 26 = arctan(2gr/A). Black dots show time intervals 
of r m t = 5. 


pose the problem of finding a preferable measurement 
basis without focusing solely on the ring up evolution for 
the resonator, and since this regime is not relevant to 
actual experiments. 


IV. QUANTUM JUMPS IN EIGENBASIS 


In this section we focus on understanding the expo¬ 
nential decay at long times in Figs. 2 and 3 for the adia¬ 
batic and frustrated regimes (b) and (c). The ensemble- 
averaged simulation suggests that after an initial state 


collapses to one of the two eigenstates { 110 ), | 01 )}, these 
eigenstates then become further mixed at a rate that de¬ 
pends on the relative magnitude of the resonator decay k 
and the qubit-qubit detuning A. As wc will soon see us¬ 
ing quantum trajectory simulations, this mixing process 
can be identified as stochastic quantum jumps between 
otherwise stabilized eigenstates. 

Treating these jumps semiclassically as telegraph noise, 
we surmise there must exist two unidirectional switching 
rates for randomly transitioning from the state 110 ) 
to |01) (—) or vice versa (+). The eigenstate population 
therefore should obey the simple ensemble-averaged rate 
equation 


^io — ^sw -^lo + (1 Pro), (19) 


where we used Pj q + P^f = 1. In particular, the solution 
of this equation starting with Pjg( 0 ) = 1 is 




r+ 

sw 

r^w + r s - 


r^w + r s 


= -(r s + w +r B - w ) t 


( 20 ) 


it eventually saturates at the population r+ v /(r+ v +T/ w ), 
and has an initial decay slope of T^. If the switching 
rates are equal, the solution will eventually reach the 
maximally mixed eigenpopulation of 1/2 (i.e., the center 
of the Bloch sphere in Fig. 3). We derive the switching 
rates for this model in the next section, after which we 
will describe how to simulate the quantum trajectories 
that show this switching behavior explicitly. 


A. Switching rate 

In order to calculate the rates of jumps between 
the eigenstates 110 ) and 101 ) in the slow dephasing regime 
r m <C A of Fig. 2(b,c), we will take literally the inter¬ 
pretation of the ensemble dephasing r m in Eq. (17) as 
being due to a fluctuating number of photons in the res¬ 
onator, causing a fluctuating ac Stark shift. Moreover, 
we will treat the fluctuating photon numbers n±(t) = 
n± + 5n±(t) as classical variables, with the means n± 
given in Eq. (13) and fluctuations 5n±(t) having tempo¬ 
ral correlations [1, 36] 

(6n±(t) $n±(0)> = n± e~^ /2 . (21) 

Here the upper sign corresponds to the state 110), for 
which the main qubit is practically in the state 11 ), while 
the lower sign is for |01). We will be mostly interested in 
the switching rate for the process 110 ) —> | 01 ), which 
is caused by fluctuations of n + (t); however, for complete¬ 
ness we calculate both switching rates (the upper sign in 
all equations below is sufficient to find T/^). Note that 
the decay rate of k/2 in Eq. (21) is consistent with the de¬ 
cay of classical energy fluctuations in a pumped resonator 
(in contrast to the energy decay k in an unpumped res¬ 
onator). Also note that here we neglected the oscillations 
of the correlator with frequency A r ± \ (discussed later). 














The fluctuating number of photons n± (t) causes a fluc¬ 
tuating ac Stark shift 2x n ±(t) [as follows from the dis¬ 
persive coupling of Eq. (4)], which leads to a fluctuating 
qubit-qubit detuning A + 8 A, with 8A(t) = 2\5n±(t). 
This in turn produces a fluctuating effective coupling 
g{t) between the stationary eigenstates 110) and |01), 
since they are no longer true eigenstates for the detun¬ 
ing A + (5A. The fluctuations 8n + (t) for the state 110) 
produce the coupling 

g + (t) = <0l|Aff|T0) = ~^8A (t) = -2^ X $n+(t), (22) 


while for the state |01) the fluctuations 8n_(t) are some¬ 
what different, producing 

~ 9 -(t) = <T0|<5tf|0T) = -2^ X Sn-{t), (23) 

where in the single-excitation subspace SH = 
{5A/2) (|10)(10| - |01)(01|) and 


n = A^l + (2<?/A) 2 « A (24) 


is energy difference between 110) and 101) (we omit the 
subscripts in A±, fi± and SA± for brevity). The deriva¬ 
tion of Eq. (22) is very simple when g -C A. Then 
the true eigenstates should correspond to the rotation 
angle 6 ss g/(A + 8 A) from the bare basis instead of 
the angle 9 ss g/A for 110) and |01). The additional 
angle, 88 ss —g8 A/A 2 , is the rotation g/A between 
the true and stationary eigenbases. Thus we obtain 
g = —(g/A)8A, which is Eq. (22) with fl ss A. In the 
exact derivation we can use 9 = arctan(2g/A)/2, then 
the derivative is d9/dA = —g/S} 2 , which should be equal 
to (g/U)/8A; this gives Eq. (22). 

The fluctuating effective coupling g between the eigen¬ 
states 110) and |01) leads to a gradual mixing between 
them, which corresponds to random jumps between the 
eigenstates in the approach of quantum trajectories. We 
can find the rate T sw of these jumps by starting with one 
of the eigenstates and equating T sw t with the population 
of the other eigenstate, which follows from the solution 
of the Schrodinger equation with the coupling g. Thus, 
to lowest order in g we find the switching rate 


T T = 

sw 


f g±(t') e ±in(t_t,) dt' 

J 0 


(25) 


where the brackets mean averaging over the random 
realizations of g{t). This equation formally depends 
on time t; however, there is actually no time depen¬ 
dence for sufficiently long t, for which the evolution 
can be physically described by a switching rate. This 
can be seen by expressing the square of the windowed 
Fourier transform in Eq. (25) via the (two-sided) spec¬ 
tral density S s± (u) of g±(t): T^ = f// 0 ^ Sg ± (±fi + 
w)[l — cos(cut)](Trtuj 2 )~ 1 duj. Therefore, at sufficiently long 
times Tf w = S§ ± (±fl), which does not depend on time. 
Because of the linear relations (22) and (23) between 


g±(t) and 8n±(t), their spectral densities are related as 
Sg ± ( fi) = (2xg/tt) 2 Ss n± (.ty, therefore 


r? w = (2 X g/n) 2 Ssn ± (±n). (26) 


[Note that for classical fluctuations 8n±(t) the spectral 
density is symmetric, Ss n± (— II) = Ss n ±(^)’, however, we 
keep the sign of SI in Eq. (26) to discuss the asymmetric 
case later.] The (two-sided) spectral density Ss n± can be 
found via the Wiener-Khinchin theorem [37] S$ n± (SI) = 
J^° 00 (^ n ±{t) 8n±(0)) e~ lQt dt, so that using Eq. (21) we 
obtain the switching rate 


= 2 g 2 8 x 2 n± k 2 
sw Q 2 k k 2 +40 2 ' 


(27) 


This result obviously assumes T^, -C |fl| and is not ap¬ 
plicable during the initial transient evolution due to col¬ 
lapse. 

Note that the term 8 x 2 n±/n in Eq. (27) is similar 
to the measurement-induced dephasing T m given by Eq. 
(17), but it depends on n + for T^ (or on h_ for T+,) 
rather than the combination fi + fi_/fi max in Eq. (17). In 
the case when n. j_ ~ ?1_ ss n max (which occurs when 
|A r ± xl < «) w e obtain T+, « T^ = T sw with 


ff 2 K 2 
A 2 k 2 + 4A 2 ’ 


(28) 


where we also used fl « A since g « A. Note that in 
the regime k > A [as in Fig. 2(c)] the last factor in Eq. 
(28) is close to 1, and the switching rate is rather large, 
T sw ss 2T m (g/ A) 2 , while in the regime k <C A [as in Fig. 
2 (b)] the switching rate is additionally suppressed by the 
factor (k/2 A) 2 . 

We have numerically verified Eq. (27) for the switch¬ 
ing rates by comparing the telegraph noise solution of 
Eq. (20) to the ensemble-averaged population decay ob¬ 
tained from simulating the full master equation (9) for a 
range of k/A. The comparison is plotted in Fig. 4, show¬ 
ing excellent agreement. A typical mismatch between the 
analytics and numerics is on the order of 0.1%, which is 
comparable to the inaccuracy from the numerical fitting 
procedure. Most importantly, in Fig. 4 we see a strong 
suppression of the switching rate at k/A <C 1. 

Solution of Eq. (19) with the switching rates given by 
Eq. (27) is sufficient to describe the ensemble-averaged 
evolution when the initial state is an eigenstate. If this 
is not the case, we need to include collapse of the initial 
state into the eigenbasis. In particular, for the bare ini¬ 
tial state 110), the evolution in Eq. (19) effectively starts 
with Pj g(0) = cos 2 9. As an example, the yellow line in 
the left panel of Fig. 2(c) shows such evolution, which is 
then converted back into the bare basis. While this sim¬ 
ple approach does not describe the transient (collapse) 
dynamics, it accurately describes the evolution after that. 

Our derivation for the switching rates T^, in this 
section has been based on treating fluctuations Sn±(t) 
as classical fluctuations. The quantum nature of these 
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effectively contribute a term (l/T)) to only the down 
switching rate T/j^. However, while such energy decay 
may be qualitatively similar in its effect on the excited 
population Pjq, it is intrinsically different from the eigen¬ 
state switching behavior derived here since it transfers 
the excitation to the ground state | 00 ) outside the single¬ 
excitation subspace, instead of switching to | 01 ). 

We also note that the mechanism discussed here of 
switching between the eigenstates 110 ) and 101 ) is phys¬ 
ically similar to the mechanism of “dressed dephasing” 
[38, 39], in which the role of the two-qubit coupling is 
played by the Jaynes-Cummings coupling between the 
qubit and resonator. 


FIG. 4. (Color online) Switching rate vs. cavity decay 
rate. Blue solid line: switching rate analytics r sw /r m = 
2(g/n) 2 /[l + (20/k) 2 ] as a function of the ratio k/A. Red 
boxes: numerical switching rate obtained from solving the 
master equation and extracting the decay rate r aw from fit¬ 
ting Pjo(t) [as shown in the right panels of Figs. 2(b,c)] with 
Eq. (20), assuming F^, = A typical relative difference 
between the numerical and analytical results is about 10~ 3 , 
which is comparable to an inaccuracy from the fitting pro¬ 
cedure. The used parameters are: g/A = 1/10 and 1/20, 
n = 10, A r = 0, and x/A = 10 -3 . 


fluctuations leads to an asymmetric spectral density 
[36] S Sn± (Sl) = n±K/[(K/2) 2 + (Q - A r =f y) 2 ]. In¬ 
serting this formula into Eq. (26), we obtain T^ = 
2{g/VL) 2 {8x 2 n±/ n)n 2 /[n 2 + 4(±fl — A r qp ,\) 2 ], which in¬ 
troduces a slight correction compared to Eq. (27). Phys¬ 
ically, this formula says that if the extra photon energy 
—A r q=y is positive, this helps the switching process with 
increase of energy, and vice versa. Even though this cor¬ 
rection is very minor in the typical case, our numerical 
results using the master equation confirm the presence of 
this correction. However, our numerical results are more 
consistent with the combination n 2 + 4(±fi — A r ± y ) 2 
in the denominator of the equation. This combination 
means that the process depends on the extra photon en¬ 
ergy —A r ± y after the switching (which changes the res¬ 
onator frequency by ±2y) instead of the extra energy 
—A r =p y before the switching. Note that the logic of our 
derivation cannot correctly take into account the change 
of the resonator frequency during switching. Also, nu¬ 
merical results for some parameters are not consistent 
with the combination — A r ±y as well (indicating a possi¬ 
ble presence of a parameter-dependent coefficient in front 
of ±y). Therefore, we are confident only in the correction 
of Eq. (27) due to A r , 


= 2 g 2 8 y 2 n ± k 2 

sw Q 2 k k 2 + 4(±fl — A r ) 2 ’ 


(29) 


omiting the dependence on y in the denominator. 

Note that if additional environmental dephasing T e is 
included in the master equation, it will contribute a simi¬ 
lar term of ( 2 T e )(g/S2) 2 to both the up and down switch¬ 
ing rates equally. Environmental energy-decay will also 


B. Quantum trajectory simulations 

In order to justify our understanding of the exponen¬ 
tial decay in Fig. 2(b,c) as resulting from quantum jumps, 
we must go beyond the master equation in Eq. (9) and 
consider more detailed quantum trajectories [30, 40-43] 
(which have been confirmed experimentally with super¬ 
conducting qubits [44-47]). In this approach we simu¬ 
late individual realizations of the evolution due to mea¬ 
surement, rather than ensemble-averaged dynamics. In 
particular, in this case there is no measurement-induced 
dephasing (i.e., a change of the qubit phase); instead, the 
gradually-acquired information obtained from measure¬ 
ment causes continuous stochastic “attraction” to the 
states | 0 ) and | 1 ) of the measured qubit (random motion 
along the meridians on the Bloch sphere). After ensemble 
averaging, these two evolutions produce the same effect, 
but in each individual measurement the effects are dras¬ 
tically different. Most importantly, using the approach of 
trajectories we simulate actual experimental realizations, 
which is impossible using the master equation. 

Since the full quantum trajectory simulation [30, 41] 
of our system is very difficult computationally, we per¬ 
formed the simulation only in the regimes of Figs. 2(a) 
and 2 (c), i.e., assuming the “bad cavity limit”, n 
(A,T). In this case the full simulation can be replaced 
with the simple quantum Bayesian approach [31, 42]. For 
further simplification we assumed that the resonator is 
driven practically on resonance, |A r ± y| -C k, and the 
setup in Fig. 1 uses a phase-sensitive amplifier, which am¬ 
plifies and outputs the optimal quadrature I(t), carrying 
information about the qubit state (the use of a phase¬ 
preserving amplifier can be described by introducing a 
limited quantum efficiency, rj < 1 / 2 ). 

The simulations have been performed in the standard 
quantum Bayesian way [ 6 , 31, 42, 44, 48], restricted to 
the two-qubit single-excitation subspace, i.e., we simu¬ 
late evolution of the density matrix with elements pio.icn 
/?oi,oii and £> 10 , 01 , using the bare basis. In brief, at each 
(small) time step d t, the unitary evolution due to the 
two-qubit Hamiltonian H q + H qq [see Eqs. (2) and ( 6 )] is 
interleaved with the evolution due to measurement, cal¬ 
culated in the following way. First, the value of the out- 
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FIG. 5. (Color online) An example of quantum jump (switch¬ 
ing event) between eigenstates 110) and 101), obtained in the 
quantum trajectory simulation for g/A = 1/20, F m /A = 
2.5 x 10~ 4 and k A (bad cavity regime). The eigenstate 
Bloch coordinate z e = Pjq — Pgy noisily hovers near ±1, ex¬ 
cept when it rapidly jumps between the eigenstates on the 
timescale of r^ 1 . Averaging these random jumps produces 
the decay observed in Fig. 2(b,c). The physical picture of 
these jumps is used later to calculate the measurement error. 

put signal I(t ) (averaged over the duration df) is picked 
randomly from the probability distribution 

e -(I-h) 2 /2D e -(I-I 0 ) 2 /2D 

p{I) =Pio,io(t) {2nD)1/2 + Pot,or W { 2 nD) i/2 ’ 

(30) 

where Ji = 1 and Jo = — 1 correspond to the bare qubit 
states 110) and |01), and the variance of the Gaussians is 
D = r/dt with the distinguishability time r = (2?7r m ) -1 
defined in Eq. (18). After picking a random value of J, 
the density matrix is updated using the relations 

Pio,io(t + df) = pio,io(t) exp[-(J-Ji)df/2r] 

Poi,oi(t + d t) poi,oi{t) exp[—(J - J 0 ) dt/2r]’ 

Pio,oi(^ + dt) _ Pio,oi(t) e ??r '"- >d4 

\J Pw,w(t + df)poi,oi(t + d t) yj pio,io{t)poi,oi{t) 

where pio,io + Poi.oi = 1 and T = r m + T e may include 
additional environmental dephasing T e . For clarity, in 
what follows we assume 77 = 1 and F = F m . For a suf¬ 
ficiently small d t, this random sampling and state up¬ 
date procedure approximates continuous stochastic tra¬ 
jectories for the two-qubit state p(t), as well as for the 
normalized readout 1 (f) = z(t) + £(f) that tracks the 
bare population difference z(t) = pio.ioW — Poi.oiWi up 
to additive white noise £(f) with a constant (two-sided) 
spectral density S = r. While the simulations are per¬ 
formed in the bare basis, the resulting density matrix p 
can be easily converted into the eigenbasis. In particular, 
we are interested in tracking the eigenbasis populations 
Pjq = pj g ig and Pqi = pgj qj besides the bare basis pop¬ 
ulations P 10 = pio.io and P m = p 0 i,oi- 

With these simulations, we can analyze transition be¬ 
tween the regimes of Fig. 2(a) and 2(c), discussed in Sec. 
III. In the textbook regime (a), r m A, we observe that 


an initial state (in the single-excitation subspace) gradu¬ 
ally collapses to either the bare state 110 ) or | 01 ) within 
the time scale of (few times) T” 1 , with rare transitions 
between the bare states at long time (note that 77 = 1 and 
ft >> A,r m ). In contrast, in the regime (c), r m -C A, we 
find from the simulations that the individual trajectories 
indeed collapse to the eigenstates 110 ) and | 01 ) at the 
same time scale T” 1 , as expected from the master equa¬ 
tion simulations. We also observe the expected random 
quantum jumps between these eigenstates at longer time 
scales. An example of such a quantum jump obtained 
from the simulations is presented in Fig. 5, showing the 
eigenpopulation difference z e = Pjo — Pqi switching from 
1 to —1. The typical “width” of the jump is compara¬ 
ble to r- 1 , though its central part can be significantly 
shorter. In between these random jumps the states re¬ 
main close to eigenstates (though sometimes with “at¬ 
tempts” of jumps), confirming the assumptions made in 
the telegraph noise model of the switching. 

Ensemble averaging of the jumps produces the grad¬ 
ual decay of the population shown in Fig. 2(c). We 
have checked numerically that the averaging of the quan¬ 
tum trajectory results coincides with the master equa¬ 
tion results, thus also confirming the formula (28) for the 
switching rate in the regime k A. Note that in our 
trajectory simulations so the switching can 

be characterized by a single rate r sw . 

We thus numerically confirm our intuitive understand¬ 
ing of the collapse to the eigenstates and rare switching 
between them when F m <C A. Note that in the simu¬ 
lated regime when also A <§C ft, the switching rate is rela¬ 
tively large, F sw « 2r m (g/A) 2 , as follows from Eq. (28). 
This can be understood as because the fast resonator de¬ 
cay ft allows each pump photon to probe only the bare 
states of the first qubit before escaping to be collected. 
In other words, the “incremental” measurement informa¬ 
tion is still sensitive to the bare basis, even though the 
relatively fast interqubit dynamics, A r m , causes the 
collapse to occur in the eigenbasis (this is because by 
the time “significant” information is collected, the eigen¬ 
basis emerges as more relevant physically). The tension 
between the different bases for the measurement in this 
“frustrated” regime leads to a relatively large switching 
rate. In contrast, when ft <C A (and still F m <C A), each 
photon in the resonator has sufficient time to feel the 
two-qubit dynamics averaged over the fast oscillations 
A. Therefore, even the “incremental” information in the 
measurement is sensitive to the eigenbasis, thus making 
it very stable and correspondingly reducing the switching 
rate r sw . This is a qualitative physical interpretation of 
the reduction factor ft 2 /(ft 2 + 4A 2 ) in Eq. (28) from the 
point of view of quantum trajectories. 

Note that this interpretation is very different from the 
physical picture used in our derivation of T sw in Sec. 
IV A, in which the reduction factor ft 2 /(ft 2 + 4A 2 ) came 
from non-zero correlation time of the fluctuations 8n(t). 
Actually, that picture was based on “fake” trajectories 
for Sn(t) and was not capable of producing collapse and 
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switching. However, it was capable of describing the 
ensemble-averaged dynamics, from which we derived T sw 
indirectly, by associating the ensemble-averaged dynam¬ 
ics with the physically correct picture of quantum jumps. 
The difference between the two pictures is that quan¬ 
tum trajectories in this section describe actual homodyne 
measurement, while in Sec. IV A we implicitly assumed 
a power (photon number) measurement right after the 
resonator. The two pictures produce the same ensemble- 
averaged dynamics because of the causality principle, but 
describe very different evolutions in individual realiza¬ 
tions of the measurement. 


Since the causality principle is not entirely trivial, let 
us discuss it in a little more detail. Classical causality re¬ 
quires that an experimenter’s action at the present time 
cannot affect anything in the past. More specifically, the 
choice of a particular action cannot affect the past. For 
example, such a choice cannot affect the evolution of an 
object that has interacted in the past with another ob¬ 
ject, which is now accessible to the experimenter. How¬ 
ever, as we know from the Einstein-Podolsky-Rosen-Bell 
paradox, this classical causality principle does not work 
in quantum mechanics, leading in particular to subtle 
“delayed choice” experiments. As a recent example, for 
a qubit continuously measured in a circuit QED setup, 
the choice of a measured microwave quadrature (selected 
with a phase-sensitive amplifier) can dictate the qubit 
evolution either along meridians or along parallels of the 
Bloch sphere [31], even though this choice affects the mi¬ 
crowave only after its interaction with the qubit (this pre¬ 
diction has been confirmed experimentally [44]). Thus, 
an experimenter’s choice in the present may affect the 
past. However, such passing of information into the past 
cannot be “useful”, in the sense that another experi¬ 
menter in the past cannot extract information about the 
later choice (otherwise it would be possible to send classi¬ 
cal information to yourself in the past; this requirement is 
often called “no signaling”). Technically, this limitation 
is caused by necessarily random results of the measure¬ 
ment: randomness saves causality. In the above example, 
we can force the qubit retroactively to move along merid¬ 
ians or along parallels, but we cannot control whether the 
qubit will move right or left (up or down), making it im¬ 
possible to distinguish the two cases (without using addi¬ 
tional information about the microwave measurement re¬ 
sult). Because of the no signaling requirement, the choice 
of the measurement cannot affect the ensemble-averaged 
evolution in the past (i.e., averaged over the random 
measurement result in the present), because otherwise it 
would be possible to extract classical information about 
the choice. (A similar argument leads to the no-cloning 
theorem [49, 50].) Thus, the causality principle in quan¬ 
tum mechanics does not forbid an experimenter to affect 
evolution in the past; however, the ensemble-averaged 
evolution in the past (averaged over randomness) cannot 
be affected by an experimenter’s choice (see also [31]). 


V. QUBIT MEASUREMENT ERROR 

The switching events (quantum jumps) contribute to 
the measurement error, which we discuss in this section. 
Here we consider only the realistic case T m <C A, when 
the eigenbasis is preferred for logical encoding over the 
bare basis. The goal of this section is to find the min¬ 
imum error, determined by the switching rate T^, (the 
analysis is very similar to the error limited by the energy 
relaxation time Tj). For simplicity we do not consider 
transients, assuming that the measurement occurs in the 
steady state. 

We assume that for the readout the information¬ 
carrying quadrature I(t) of the output signal is integrated 
over the measurement duration t, producing the averaged 
output 

/(*) = IjfVx (32) 

and then this value is compared with the threshold J t h to 
produce a binary readout of “0” or “1.” (More advanced 
signal processing of I (t) can moderately improve the mea¬ 
surement fidelity [22, 51]; we consider the straightforward 
integration (32) for simplicity.) If our goal is to distin¬ 
guish the states 110 ) and 100 ), then the probabilities of 
misidentifying these states are 

/ Ith _ _ /*00 

P(/|T 0 )d/, P^) = P(7|00) d/, (33) 

-oo J /th 

respectively, where P(71 10) is the probability density of 
obtaining the result 7 when the initial state is 110 ) and 
P(7 | 00) is the analogous probability for the initial state 
100). The total measurement error is the average of the 
two errors, 

Perr=\\P^+P^}. (34) 

Figure 6 shows example histograms for P(7 | 00) and 
P(J | 10), obtained by simulating 100,000 quantum tra¬ 
jectories, as discussed in the previous section for t/r = 7, 
g/A = 1/10, r m /A = 10~ 3 , and r] = 1 (in this case 
r sw r = 9.6 x 10~ 3 ). As in the previous section, we use 
the normalization in which an ideal single-qubit measure¬ 
ment corresponds to I = I\ = 1 for the state | 1 ) and 
7 = Jo = — 1 for the state |0). As seen from Fig. 6, the 
probability distribution P(71 00) is a Gaussian centered 
at 7 = — 1 , while P(71 10) has a significant “tail” (red 
shaded region), caused by switching events. Also, the 
Gaussian part of P(7110) is centered at a value slightly 
smaller than 1 (this shift is practically not visible) be¬ 
cause the eigenstate 110 ) has a small contribution from 
101). The shape of the histograms is discussed in the 
appendix. 

The errors PirV, Pei ), and P err depend on the choice of 
the threshold It h- Obviously, the total error P err is min¬ 
imized when the threshold is set such that P(/th 1 10 ) = 


12 


P(I) 



FIG. 6. (Color online) Signal histograms P(I \ 00) and 
P(I | 10) for the integrated quadrature 7, given the initial log¬ 
ical eigenstates 100} and jlO). Shown are the binned readouts 
for 100, 000 trajectories simulated as in Fig. 5 for duration 
t/r = 7, with g/A = 1/10, F m /A = 10~ 3 , and g = 1. Note 
that P(I | 00) is a Gaussian centered at Jo = —1, but P(I | 10) 
is a slightly shifted Gaussian [centered at 1 — 2 (g/A) 2 instead 
of Ji = 1], with a significant extended “tail” (red shaded re¬ 
gion) caused by quantum jumps. The overlap of the two Gaus- 
sians decreases with integration time t, but the histogram 
overlap due to the tail increases with t, thus preventing per¬ 
fect discrimination. 



FIG. 7. (Color online) Simple analytics for the measure¬ 
ment error P e rr as a function of an integration duration 
t, normalized by the distinguishability time r. Orange 
dashed line shows the monotonically decreasing separation 
error from integrated white noise. Green dashed line shows 
the linearly increasing error F.^t/4 from switching events, 
for which we choose = 10~ 3 . The blue solid line 

shows the total measurement error, which has a minimum of 
P e ™ m « (F^r/2) ln(0.6/r a - w r) at the optimum time t op t/r « 
2 ln(l/4F s _ w r). 


P (/th | 00). However, in most of this section we will as¬ 
sume the symmetric threshold, Ith = 0. This is done 
for simplicity and also because, as we will see later, the 
use of the optimal threshold decreases the error insignif¬ 
icantly (with a typical relative improvement of < 3%). 
Also note that we will discuss the error for distinguish¬ 
ing the eigenstates 110) and |00) as optimal for logical 
encoding; the corresponding error for distinguishing the 
bare states 110) and |00) has an additional contribution, 


cos20)/2].) For a small g/A this correction is small, and 
we will neglect it below. 

The separation error rapidly becomes very small: ICG 2 
for t = 5.4 r, 1CP 3 for t = 9.5 r, and 10~ 4 for t = 13.8 r. 
However, the switching process 110) —l |01), occurring 
with the rate T^,, adds a contribution to the error Pirr 
that increases in time nearly linearly, 



P £,bare « H)) + (ff/A) 2 , (35) 

because of the initial collapse of the bare state 110) into 
either the eigenstate 110) or |01). 


A. Error contributions 


In the absence of switching, T^, = 0, the error steadily 
decreases with integration time because the variance 
cr 2 = r/t of the Gaussians in Fig. 6 decreases with time 
t [the distinguishability time r is given by Eq. (18)]. In 
particular, for Ah = 0 this “separation” error is 


p(o) 

err,sep 

pi 1 ) 

err,sep 


1 - erf(v / t/2r) , 

1 — erf(cos(20)-\/i/2T) , 


(36) 

(37) 


where cos(20) « 1 — 2(g/A) 2 comes from the difference 
between the eigenbasis and the bare basis. (For the op¬ 
timal threshold both errors will contain erf[y7/2r(l + 


so that the total error becomes 


P>, 


1 — erf( v / t/2r) 


-r 

4 sw 


t, 


(39) 


where we used cos(20) « 1, J th = 0, and (T^, + F^f <C 
1. More accurate calculations (in particular, taking into 
account double-switching trajectories, proper convolu¬ 
tion of switching and noise, and effects of 6) are presented 
in the appendix. Note that it is easy to understand the 
factor 1/2 in Eq. (38) by saying that the initial state 
110 ) will be misidentified only if the switching event oc¬ 
curs before the middle of the integration time, so that 
the erroneous state is integrated for a longer time than 
the correct state. (A better interpretation of this factor 
via symmetry of the convolution is discussed in the ap¬ 
pendix.) We also note that Eqs. (38) and (39) can also 
describe the error for a single-qubit measurement that ac¬ 
counts for energy relaxation, with replaced by T) -1 . 

In Fig. 7 we illustrate the decomposition of the eigen¬ 
basis measurement error P err of Eq. (39) into its two parts 
for r- w r = 10~ 3 . The orange dashed line shows the error 
contribution from the integrated white noise (separation 
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noticeably shifts the curve down at times after the mini¬ 
mum. [Most of the difference between the green and blue 
lines at the times after the minimum can be obtained by 
simply replacing T~ w t in Eq. (38) with 1 — exp(—r^t); 
the remaining difference is mainly due to double switch¬ 
ing events.] As we see, the analytics of the appendix is 
closer to the numerical results than the simple analytics, 
but the difference is minor. This difference becomes even 
smaller for smaller ratios r^/r m . 


B. Error minimized over time 


FIG. 8. (Color online) Simulated measurement error PerV for 
misidentifying the initial state 110), using the discrimination 
threshold of 7th = 0. Red solid line: measurement error ob¬ 
tained by binning the integrated readouts for M = 1, 500, 000 
individual quantum trajectories, with g/A = 1/20, A/r m = 
2000, and g = 1, in the bad cavity regime (k A). Error 
bars show the standard deviation of [Pjrr (1 — PjrV )/M] 1 ^ 2 . 
Green dashed line: simple analytics that includes only single 
quantum jumps, Eq. (38). Blue dot-dashed line: refined an¬ 
alytics that includes single and double quantum jump events 
(see the appendix). 


error), which monotonically decreases with integration 
time. The green dashed line shows the error from switch¬ 
ing events, which linearly increases with integration time. 
The combination of these two opposing effects in the total 
measurement error (blue solid line) produces a minimum 
error P“ r ln at an optimum time f opt , which we discuss in 
the following section. 

To verify that this simple approach adequately mod¬ 
els the measurement error, Fig. 8 shows a comparison 
of the measurement error Pi r V for the initial state 110) 
calculated in three ways: using quantum trajectories, us¬ 
ing the simplified description (38), and using the more 
accurate analytics discussed in the appendix. The quan¬ 
tum trajectory method has used M = 1,500,000 indi¬ 
vidual trajectories initialized in the eigenstate 110) for 
g/A = 1/20, A/T m = 2000, and rj = 1. Each trajectory 
consists of 2x 10 5 time steps of size dt/r = 10 -4 . For each 
trajectory we calculate 7(f) via Eq. (32) and compare it 
with the threshold 7 t h = 0. The error P e rr is then the 
fraction of trajectories with I < 7 t h, which is shown by 
the solid red line in Fig. 8. Error bars show the standard 
deviation [Pi r V (1 — PirV )/M] 1 / 2 for a few representative 
points. For comparison, the dashed green line shows the 
simple analytics (38) with = 2T m (g/VL) 2 [see Eqs. 
(27) and (28)], so that F^r = 1/404. The dot-dashed 
blue line shows the more refined analytics described in 
the appendix that include double switching events, as 
well as the proper offset of the Gaussian by cos(20). This 
offset slightly shifts the curve up at times before and near 
the minimum, while the more accurate account of jumps 


The simplified model (39) for the measurement er¬ 
ror has a first contribution that is rapidly decreasing in 
time, and a second contribution that is slowly increas¬ 
ing in time (Fig. 7). Therefore, it has a minimum that 
is reached at an optimal time 7 opt . The minimum error 
P e ™ n at this optimal time should be determined by the 
product F s - W r, since this is the only dimensionless param¬ 
eter in the model. The optimal measurement duration 
can be found via the equation dP err /df = 0, whose so¬ 
lution is a product-log function, which has the recursive 
form 


7 op t = t In 


2 / 7T 


(T sw T) 2 (t opt /T) 


The corresponding minimum error is 


(40) 


r- 7 

pmin ^ sw' 
err ~ 2 


Tsw^opt 


r- T n 

In 0.6, 

2 r sw r 


(41) 

where C ss e-\/2/7r-y/r/7 0 pt actually depends on r sw r, 
but sufficiently weakly: C € [0.43, 0.74] for rF^ £ 
[10 -6 , 10“ 2 ]. In deriving the first relation in Eq. (41) 
we used the approximation 


1 — erf (a;) 
2 


exp(— x 2 ) 
2y/ir x ’ 


x » 1 . 


(42) 


Note that in the usual case tT^ -C 1, the main con¬ 
tribution to P™ r ln comes from the second term in Eq. 
(39), with the relative contribution from the first term 
being 2r/7 op t. The ratio t op t/T can be estimated as 
topt/r « 2 ln(C/eF s - w r) « 2 ln(l/4r^r). 

Note that at very long times, t > (T^ + T+,) -1 , 
the simplified model in Eq. (39) becomes inapplicable. 
Moreover, the reverse switching events with the rate 
r+y will eventually produce the integrated output signal 
I -t (T+, - r- w )/(r+, + r^,) for the initial state |10), 
which thus can be distinguished with certainty from the 
initial state 100) (a similar situation was discussed for the 
measurement of phase qubits in Ref. [21]). However, such 
long integration times, t 3> l/(r/ w +r/ w ), are impractical 
even if we assume the absence of the energy relaxation, 
so we do not consider P e „(t) for these long times. 
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C. Optimized threshold 


Let us augment the simplified model (39) by introduc¬ 
ing an arbitrary threshold I t h ; then the error becomes 

p _ l-erf[(l +J th )v / t/2r] 

err ~ 4 

+ l-erf[(l-/,„)yi7^] +i +fa, r _^ (43) 

Choosing a slightly negative J t h decreases the error be¬ 
cause of the contribution from the last term (this is also 
obvious from Fig. 6 since the histogram for the initial 
state 110) has a long tail). The optimal threshold 1°^ 
and optimal time f opt can now be found from the system 
of equations, dP err /d/th = 0 and dP crr /dt = 0. These 
equations are rather lengthy, but in the case t op t/r 3 > 1 
lead to a simple relation exp(— Ithtopt/r) = v/3- There¬ 
fore, the optimal threshold is only slightly different from 
zero 


rOpt 

J tii 


In 3 r 

2 t Q pt 


-0.55 


^opt 


(44) 


while the optimal time f 0 pt does not change significantly 
compared with Eq. (40). 

The optimal threshold 1°^ can also be obtained in the 
following crude way. Using Eq. (43), let us calculate the 
first and second derivatives of P err over J t h at the point 
Ith = 0. This is simple and gives dP orr /d/ t h = T~ w t/A, 
d 2 P err /d/ t 2 h = (27r) _1 / 2 (t/r) 3 / 2 exp(— t/2r). Then, as¬ 
suming a parabolic dependence P err (7th), we find the op¬ 
timal threshold as 1°^ = — (dP err /d/ th )/(d 2 P err /d/ 2 h ), 
which, also using Eq. (40), gives 1°^ = — T/2t opt . There¬ 
fore, this crude derivation does not reproduce the re¬ 
sult (44) exactly, but is still quite accurate. The er¬ 
ror decrease due to optimization of Ith can then be 
found from the same parabolic approximation as <5P™ n = 
(dPerr/d/thKt^V 2 ! which gives 

P“ in (J th = J t ° h pt ) - P“( n (/ th = 0) « (45) 

Since we do not expect a significant change of P err due 
to a slight shift of t opt in this double-optimization proce¬ 
dure, we can simply replace the term (1/2) in Eq. 
(41) with (7/16) T^t. Therefore, in this crude derivation 
the error P™ r ln optimized over both time and threshold 
is still given by Eq. (41) with a modified value of C, 


Copt = e" 1/8 C « 0.88 C ~ 0.5. (46) 

The relative decrease of P e “ ln due to the threshold op¬ 
timization is approximately [8lnfC/T^T)] - , which is 
about 3% for T^r = 10 -2 and smaller for smaller values 
ofr s - w r. 

By solving the optimization problem numerically over 
a wide parameter range £ [10 _6 ,10~ 2 ], we have 

confirmed that the threshold optimization changes C [de¬ 
fined via Eq. (41)] by a nearly constant factor C op t/C ~ 


0 .88, producing the range C op t £ [0.37, 0.65]. Corre¬ 
spondingly, this produces a nearly insignificant relative 
correction of [1.0%, 3.2%] in the minimum error -P™ 11 over 
this parameter range. The denominator 16 in Eq. (45) in 
numerical results is found to be close to 15. Similarly, we 
confirmed that Eq. (44) is satisfied quite well: instead 
of the factor 0.55, we numerically find 0.55-0.58. This 
leads to the numerically optimal threshold I t h varying 
only within [—0.023, —0.081] over this same parameter 
range. 

Besides using Eq. (43) for the numerical optimization, 
we also used a modified equation, in which the second 
term is multiplied by exp(—T^t), and in the third term 
T~ w t is replaced with 1 — exp(—T^f). This practically 
does not change the above mentioned results, except 
that it slightly lowers C: for the same parameter range 
T- W r £ [10 -6 ,10~ 2 ] it is C £ [0.43, 0.64] and correspond¬ 
ingly C 0 pt € [0.37, 0.57]. Note that the factor 4 in the 
mentioned above approximation t opt /T ~ 2 ln^/dr^r) 
varies within the range [6.1, 3.1] for zero threshold and 
the same parameter range, and within [5.3, 2.7] for the 
optimal threshold, increasing f opt by about 0.3r com¬ 
pared with the zero-threshold case. 

We emphasize that the main result of our analysis of 
the threshold optimization is that the optimal threshold 
J)] 1 / 1 is close to the symmetric point Ith = 0 and that 
the benefit of the optimization in decreasing the mea¬ 
surement error is insignificant. This justifies the use of 
Ith = 0 in the analysis. This also shows that it is mean¬ 
ingful not to perform the threshold optimization in an 
experiment and instead use the symmetric point. Note 
that this conclusion also applies to the case of qubit en¬ 
ergy relaxation, which is analyzed in the same way. 


D. Measurement error due to neighboring qubit 


As discussed above, the minimized measurement error 
can be approximated as 


pmin 

err 


r T 

sw 


In 


C 


(47) 


where C ~ 0.6 for the symmetric threshold or ~ 0.5 for 
the optimal threshold, r = (2??r m ) -1 is the distinguisha- 
bility time, and 


1 

Ti 


2 g 2 8 x 2 n+ k 2 _ 

A 2 + 4 g 2 k k 2 + 4A 2 + 16<7 2 


(48) 


is the switching rate [see Eq. (27)]. Since in this paper 
we are interested in the effect of the neighboring qubit, 
let us neglect the energy relaxation rate 1 /Tj. Also, let 
us use | .<7 /A| < 1 and assume |A r ± yj < k, so that 
8 \ 2 n+/n ~ T m . In this case is given by Eq. (28) and 

the measurement error is 


pmin 

-^err 


1 /3\ 2 « 2 , ac 2 + 4A 2 ’ 

2? ? \a) k 2 +4A 2 [ 11 ac 2 ( 5 /A) 2 


( 49 ) 
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Note that this is the error for distinguishing the states 
110 ) and |00), while the error for distinguishing the bare- 
basis states 110) and |00) is larger because of the collapse 
occuring in the eigenbasis when |A| r m , 

pmin ~ (JL\ 2 _l_ pmin (erp 

x err,bare ~ 2 V ^ J err " ) 

As we see, in the “bad cavity limit”, n |A|, the 
eigenbasis error, P e “ ln ~ (l/2rj)(g/A) 2 ln(Cr)A 2 /g 2 ), is 
quite large, for example, for g/A = 1/10 and t) = 0.2 
we obtain P e “ ln ~ 6%. The bare-basis error P^bare is 
bigger by 0.5%, which is not significant because P e 7 r m is 
so big. 

This may look dangerous for the quantum processors 
based on superconducting qubits with “always-on” inter¬ 
action between the neighbors. Fortunately, typical exper¬ 
imental systems do not operate in this bad cavity limit; 
in more realistic parameter regimes the switching pro¬ 
cess is strongly suppressed and therefore the measure¬ 
ment error due to the presence of a neighboring qubit 
is relatively small. For example, for k~ 1 = 20 ns and 
A/27T = 0.5 GHz, the switching rate is approximately 
10~ 4 r m (g/A) 2 , so that for g/A = 1/10 and rj = 0.2 we 
obtain a very small measurement error from the neigh¬ 
boring qubit, P™ n « 2 x 10 -5 . However, the bare-basis 
error is still significant, P™ r “ are ~ 0.5%, which means 
that the bare basis is inappropriate for encoding the log¬ 
ical information. Since the eigenbasis is also beneficial 
for logic operations and idling [21], this makes it an un¬ 
ambiguously optimal choice for encoding quantum infor¬ 
mation. 

Note that since the switching processes are strongly 
suppressed in the regime when k <C |A|, the qubit 
measurement remains accurate even when a neighbor¬ 
ing qubit is detuned only moderately, \A/g\ ~ 3, as long 
as the eigenstates are used for encoding. This fact may 
simplify the design of quantum processors in which “fre¬ 
quency crowding” may present a problem. 

The switching between the eigenstates 110) and 101) 
can be observed experimentally. (For this purpose it is 
better to use the jumps 101) —>• 110), which can be easily 
distinguished from energy relaxation events.) For exam¬ 
ple, for K -1 = 10 ns, g/2n = 30 MHz, A/27r = 100 MHz, 
and r m /27r = 20 MHz (corresponding to |x|/27r = 2 MHz 
and n = 10), the expected switching time is about 
r-J ~ 10 ns. 


VI. CONCLUSION 

We have investigated the measurement error of a su¬ 
perconducting transmon qubit in a circuit QED setup 
caused by the coupling g to a detuned neighboring qubit 
(or a bus resonator), focusing on the effects of the corre¬ 
sponding “tail” population (g/A) 2 . When the ensemble¬ 
dephasing rate due to measurement is much faster than 
the qubit-qubit detuning, the system collapses to the bare 


energy states as one would expect for a textbook pro¬ 
jective measurement. However, in the more physically 
relevant regime with the ensemble-dephasing rate much 
slower than the detuning, the system instead collapses 
to the joint qubit-qubit eigenstates, which are also favor¬ 
able for quantum operations and idling. As such, these 
qubit-qubit eigenstates are the most appropriate states 
for high-fidelity logical encoding in realistic parameter 
regimes. 

We have shown that in regime where joint eigenstates 
are preferred, the excitation can randomly jump between 
these eigenstates while the qubit is being measured. In 
between these random jumps, the two-qubit state is prac¬ 
tically pinned to one of the eigenstates. We have derived 
the rate of the jumps by using a semiclassical model of 
fluctuating ac Stark shift. The obtained analytical result 
for the switching rate has been confirmed by comparison 
with numerical solution of the master equation, for which 
the ensemble-averaged jumps lead to a gradual decay of 
the initial eigenstate population. The random jumps pro¬ 
duce a contribution to the measurement error probability 
that increases almost linearly with integration time in a 
way qualitatively similar to the error from energy-decay 
processes. 

The switching rate for these random jumps depends 
on the relative magnitude of the resonator decay and the 
qubit-qubit detuning. For quickly decaying resonators, 
the switching rate is significant and produces the mea¬ 
surement error exceeding (g/A) 2 by several multiples. 
However, for more slowly decaying resonators, as is more 
typical experimentally, the minimized measurement er¬ 
ror becomes essentially negligible for eigenstate encoding, 
while the error for bare basis encoding is still significant 
and exceeds \(g/A) 2 . 

For the purposes of this study, we have used a static 
threshold for digitizing the continuous quadrature read¬ 
out. We note that more sophisticated discrimination 
schemes may be able to take advantage of the addi¬ 
tional information contained in the continuous readout 
to partially correct for the switching contribution to the 
measurement error. Generalizing our analysis to mul¬ 
tiple neighboring qubits with simultaneous multi-qubit 
measurement may also be interesting for future research. 
Another possible generalization is the analysis of deco¬ 
herence for multiqubit states that involve the neighbor¬ 
ing qubit, which are not supposed to be affected by the 
measurement, but are actually influenced by the switch¬ 
ing dynamics. We emphasize that the quantum jumps 
between the eigenstates predicted in this paper could be 
measured experimentally using existing superconducting 
qubit technology. 
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Appendix: Measurement error derivation 

In this appendix we provide a more complete deriva¬ 
tion of the measurement error, assuming the regime with 
r m -C A, where the two-qubit eigenstates are the op¬ 
timal logical states. Our derivation will consist of two 
parts. First, neglecting transients for the resonator and 
assuming abrupt switching events, we will calculate the 
histograms P(I | 00) and P(I | 10) for the integrated mea¬ 
surement output I = f 0 I(t') dt'/t, corresponding to the 
initial states |00) and 110), respectively (we use the word 
“histogram” instead of “probability distribution” as a 
shorter term). Second, we will impose a discrimination 
threshold I t h on these histograms to compute the prob¬ 
ability of error according to our definition in Eq. (34). 


We note that the assumption of an uncorrelated dynam¬ 
ics of £(i) and z(t) is in general not good in the quan¬ 
tum Bayesian approach. For example, for single-qubit 
Rabi oscillations, the correlation between £(t) and z(t) is 
as strong as autocorrelation for z(t) [52], However, for 
the dynamics with rare switching events this assumption 
should be sufficiently good because the correlation be¬ 
tween £(i) and z(t) is most important only in the vicin¬ 
ity of switching events, which occupy a small fraction 
of the total integration time. The approximation (A.l) 
also neglects transients of duration ~ k -1 at the start 
of the measurement and near switching events, implying 
t k -1 . 

The approximation of Eq. (A.l) permits us to calculate 
each histogram for I in a simple way as a convolution be¬ 
tween a histogram for the population difference Ztot and 
the Gaussian white noise distribution for £. If we ad¬ 
ditionally assume that in the single-excitation subspace 
the state is always pinned to an eigenstate, with abrupt 
jumps between the eigenstates (in particular, this implies 
Pol io = 0), then the histogram for z is determined by the 
histogram for the eigenstate population difference z e , 

1 /'* 

2 = cos(26>) z e , z e = - [P^t') - P m (t')]dt', (A.5) 

t J o 

via the conversion factor cos(20) with rotation angle 
2d = arctan(2g/A). We will now calculate the his¬ 
tograms corresponding to the specific initial populations 
Poo(0) = 1 or Pjq(0) = 1, which are the optimal logical 
states for discrimination. 


1. Readout histograms 

Our primary assumption for obtaining the readout his¬ 
tograms is that we can separate the integrated normal¬ 
ized measurement output I into two approximately un¬ 
correlated terms (signal and noise) 

I(t) = Zto t (t) + €(t). (A.l) 

The first term, f t ot = z + Z, is the total integrated 
bare-population difference between the ground and ex¬ 
cited states of the main qubit, which includes a part z in 
the single-excitation subspace, as well as a part Z outside 
this subspace, 


a. Ground-state histogram 

An initial ground state 100) remains |00) for an ar¬ 
bitrarily long time, and the corresponding output signal 
also does not change in time, since we assumed the steady 
state for the resonator. Therefore, an initial population 
-Poo(O) = 1 produces the stationary integrated coordinate 
Z = — 1 (with z = 0) and the stationary integrated total 
population difference f t ot = — 1, which implies a delta- 
function histogram P z [5 tot (t) | 00] = 5(f to t + !)• Con¬ 
volving this histogram with the Gaussian white noise in 
Eq. (A. 4) produces the histogram for the integrated mea¬ 
surement result, 


m = 7 f [P w {t’) - P m {t')]&t’, (A.2) 

1 Jo 

Z{t) = 7 f [Pn(t') - P 0 o(t')]dt' . (A.3) 

t Jo 

The second term of Eq. (A.l) is integrated zero-mean 
white noise, which is randomly sampled from a Gaussian 
distribution of variance r/t, 


P(I | 00) 


/ * 

(/+ Iff 

V 27TT 6XP 

2 T 


(A.6) 


which is the expected Gaussian distribution of the same 
width as the noise that is centered at the ground state 
normalized signal of / 0 = — 1. 


b. Excited-state histogram 






2 r 


(A.4) 


An initially excited eigenstate 110) will randomly jump 
to the eigenstate |01) at the rate T^,, as discussed in Sec- 
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tion IV, and may then randomly jump back to the orig¬ 
inal eigenstate at the rate T+,. We assume that these 
jumps can be treated as instantaneous compared to the 
integration time t, and that we can treat the eigenstates 
as stationary between these jumps. We also assume that 
the jumps obey Poissonian statistics and that the aver¬ 
age time between the jumps is long compared to typical 
integration times, T^t <C 1; therefore it will be sufficient 
to consider only zero, one, or two possible jumps per in¬ 
tegration duration t. The total histogram for the excited 
state will then be a weighted contribution of histograms 
with a definite number of jumps 

P(j | TO) =p 0 P (0) (/|IO) +piP (1) (/|10) +p 2 P {2) (i\T0). 

(A.7) 

We compute each of these histograms and their weights 
separately. The derivation is significantly easier for the 
case when T+, = T^ = T sw , so we will be starting the 
discussion with this case and then discussing more ap¬ 
proximate results for unequal switching rates. 

For Poissonian jump statistics with equal switching 
rates, the probability of having k jumps within the mea¬ 
surement duration t is Pk = e~ rawt (T sv ,t) k /k\. Since we 
consider only k < 2, we will use 

—r HW t (r s w^) —r sw t i 

p 0 = e sw , P 2 =— x — e sw , pi = 1 -p 0 -p 2 , 

(A.8) 

where p\ is chosen to adjust normalization because the 
largest neglected contribution, k = 3, also has an odd 
number of jumps, and such choice slightly improves the 
accuracy for the calculation of the measurement error 
(though this is not really important). 

In the case of unequal switching rates the exact for¬ 
mulas for pk are quite lengthy, po = W 1 " 8 " 4 , p± = 
r S w(e- r =w - e- r ^)/(P s+w _ r - w) , P2 = r s - w r+[ e -fo* + 
e- r -*(r+t - r~ w t - l)]/(r+ w - r" w ) 2 , so we will use the 
approximation 

po = e- r ™‘, P2 « r ^ s+w * 2 , pi = 1 - po - P2, (A.9) 

Note that the next-order approximation for p 2 is p 2 ~ 
ir s - w r+T 2 [l — (2rjj v -l-r+fo/S]. Also note that to linear 
order in t (then fully neglecting p 2 , as in the main text) 

Pi ~ r sw*- ( A - 10 ) 


which is a Gaussian distribution similar to the ground- 
state histogram, but centered at cos(20) « 1 — 2(g/A) 2 
that is slightly shifted from the otherwise expected mean 
of I\ = +1 by the coupling to the neighboring qubit. 
This slight shift was neglected in the main text as small. 

Now let us calculate the histogram with a single jump. 
If the jump occurs at time moment t\, then the signal 
of cos(20) is integrated for time t\ and the signal of 
— 008 ( 20 ) is integrated for time t — t i, resulting in the 
average 


= —-—- cos(20). (A.13) 

In the case Pj^ = P+,, the jump time t\ is equally 
likely at any time in the interval [0, i]. Then has 
a uniform histogram pfo(z to t 110) = [2cos(20)]” 1 in the 
corresponding interval [—cos(20), cos(20)], illustrated in 
Fig. 9(a). Convolving this uniform distribution with the 
Gaussian white noise in Eq. (A. 4) yields the measurement 
histogram 


P (1) (/| 10) 


. / + cos(20) I — cos(2 9) 

erf- . -erf- . — 

y/2r/t _ a /2 r/t 

4cos(2 0) 


(A.14) 


which is a smoothed box distribution [in Fig. 9(a) the 
smoothing is shown for t/r = 10]. This addition to 
the excited-state histogram is the dominant effect of the 
quantum jumps on the readout. Note that for t/r > 4 the 
smoothing appreciably affects only the edges of the rect¬ 
angular distribution Pz^(z tot 110 ), so that for cos( 20 ) — 
|J| > 2 y/r/t we can use PW(/| 10 ) ~ PS 1 "*(P | 10 ). 

In the case of unequal switching rates, 7 ^ F^,, the 
jump time t\ is no longer equally distributed within t; 
instead, it has the (normalized) probability distribution 
e -rr w t le -r^( t - tl )(r^ w - r B -)/(e- r -* - e _r ™ 4 ). Then 
the probability distribution for within the interval 
[— cos( 20 ), cos( 20 )] is given by the same formula, with t\ 
replaced by [1 + zSP/ cos( 20 )]t /2 and extra normaliza¬ 
tion factor t /[2 cos(20)]. This probability distribution for 
should then be convolved with the Gaussian noise to 
obtain P^^/jlO). The resulting formula is very long, so 
for simplicity we can use 


The excited-state histogram with zero jumps is similar 
to the ground-state histogram in Eq. (A. 6 ). The initial 
eigenpopulation remains stationary in this case, corre¬ 
sponding to a stationary eigenpopulation difference z e = 
1 , and thus a bare population difference of z = cos(2 6) 
and a histogram of 

Pz°\z tot 110) = 5foot - cos(20)]. (A.11) 


Convolving this histogram with the Gaussian noise in 
Eq. (A. 4) yields the measurement histogram 


P (0) (/| 10) 


/ t 

[I — cos( 20 )] 2 t 

V 27TT 6XP 

2 T 


(A.12) 


1 + v sw swy z 

Pi 1} foot 110) « - llTufnl -> l^totl < cos( 20 ), 

2 cos(2U) 

(A .15) 

and since the convolution of a linear function with the 
Gaussian noise affects mostly the vicinity of edges, we can 
use PW(/| 10 ) ~ Pfofol 10 ) for cos( 20 ) — |J| > 2 \frji. 
A little better approximation is to use Eq. (A. 14) with 
added term J(r+, — rjj v )t/4cos 2 (20) at |/| < cos(20). 

Note that if we also want to take into account the 
energy relaxation with the rate Tjj , then for the en¬ 
ergy relaxation event occurring at time ti we have = 
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FIG. 9. (Color online) (a) Left panel: schematic evolution 
of the eigenbasis population difference z e with a jump from 1 
to —1 at ti. Right panel: the histogram P- 1 ' 1 (z t ot| 10) for 
time-averaged bare-basis population difference in the one- 
jump scenario (solid line) and the corresponding histogram 
P^(/|10), which includes Gaussian noise (dashed line), (b) 
Similar schematic and histograms for the two-jump scenario. 
We assume t/r = 10, g/A = 1/10, and F^ = r,fo. 


[1 + cos(20)] t\/t — 1. Then using approximation of uni¬ 
formly distributed t\ (applicable for t/Tj -C 1) we ob¬ 
tain the uniform distribution for zf:^ within the interval 
[—1, cos(2d)[. Convolution with the Gaussian noise will 
then lead to a slightly asymmetric probability distribu¬ 
tion pW(I |T0). 

Now let us calculate the histogram with two jumps. If 
the first jump occurs at time moment t\ and the return 
jump occurs at £2, then the system spends duration At = 
f2 — ti in the “wrong” state 101), so that 

_ (2) = t-2At cog(20) _ ( A . 16 ) 

In the case F^ = T+, the probability distribution of 
time moments ti and t2 is uniform within the range 
0 < ti < t 2 < t, and therefore the normalized probabil¬ 
ity distribution for At is linearly decreasing, Pa* (At) = 
(2/t)[1 — (At)/t]. This produces the linearly increasing 
probability distribution for the integrated signal, 


P?\z tot 110) - 2 ( C0S 2 ( 2^ Ot > l^otl < cos(20), 

(A. 17) 

which is illustrated in Fig. 9(b). The convolution with the 
Gaussian white noise in Eq. (A. 4) yields the measurement 
histogram 


P (2) (/| 10) 


fCos(2 9) e -t(I- z ) 2 /2T 


’ — cos(2 6) yj 2 tTT J t 


cos(2 9) + z 

- Cl 7 

2 cos 2 (20) ’ 

(A.18) 


which we leave unevaluated here for brevity. As above, 


FIG. 10. (Color online) Comparison between analytical and 
numerical results for the measurement histograms P(I | 00) 
and P(I | 10). The green and blue lines show the same numer¬ 
ical results as in Fig. 6, but on an enlarged scale. The solid 
black line shows the analytical result [Eq. (A. 7)] for P(I\ 10), 
taking into account up to two jumps. The almost coinciding 
red dashed line shows the result with up to one jump. The 
black dashed line shows Eq. (A. 6) for P(/| 00). 


we can use approximation p( 2 \l 110) « Pz 2 \l 110) suf¬ 
ficiently far from the edges, cos(2 9) — |/| > 2 y/r/t [see 
Fig. 9(b)]. 

In the case when T^ 7^ T+,, the exact formulas are 
much lengthier because the probability distribution for 
the jump moments t\ and £2 [which is proportional to 
e^Tsw-rswH^-L)] is n o longer uniform. However, since 
|r+ w -T- w | t <C 1, we can still approximate it as a uniform 
distribution, and then Eqs. (A. 17) and (A.18) are still ap¬ 
proximately valid. Since the two-jump histogram brings 
a very small contribution to the total histogram (A.7), 
any crude approximation should be sufficient. Note that 
if the first jump was due to the energy relaxation event, 
then the return jump is impossible. 

Thus, the no-jump contribution to the total excited- 
state histogram (A.7) produces the main Gaussian shape, 
the single-jump contribution adds an extended nearly 
uniform tail, and the two-jump contribution produces a 
very small linearly increasing correction. Figure 10 shows 
on an enlarged scale the numerical (quantum trajectory) 
histograms presented in Fig. 6 (green and blue lines) and 
also shows the analytical results derived in this section. 
The solid black line shows the result for P(/|10) using 
Eq. (A.7) taking into account up to two jumps, with 
r s - w = r+, = 2r m g 2 /(A 2 + 4ff 2 ) for parameters of Fig. 6 
(so that F^r = 9.6 x 10 -3 and t/r = 7). The dashed 
red line shows a similar result using a simplified one-jump 
approach, in which p 0 = e~ r ™ t and p\ = 1 — Pq. We see 
that the difference between the results for the one-jump 
and two-jump approaches is very small, but the two-jump 
approach still agrees slightly better with the numerical 
results (blue line) for the tail of the distribution. It is 
interesting to note that the tail looks almost linearly in- 
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FIG. 11. (Color online) Time evolution of the integrated sig¬ 
nal histogram for an initially excited eigenstate 110), shown 
on a semi-log scale, for times t/r = 5, 10, 15, 20, assuming 
Fs*wT = 10 -3 . The histogram width due to Gaussian noise 
decreases with increasing integration time, analogously to the 
ground state histogram (which is not shown). Also, a tail due 
to quantum jumps between the eigenstates appears at the left 
side of the histogram; this tail grows in amplitude and flattens 
with increasing integration time. 


creasing, in contrast to the horizontal shape expected 
from the uniform distribution of P/ 1 (Aot|10). This is 
because near its edge, 7 = — cos(20) ss —1, the Gaussian 
averaging plays the major role [see Fig. 9(a)]. The black 
dashed line shows Eq. (A.6) for P(7|00); its agreement 
with the numerical results (green line) is rather trivial 
because in this case only noise is simulated numerically. 

The evolution of the tail of the excited-state distribu¬ 
tion P(7110) is illustrated in Fig. 11 (we do not show 
the ground-state histogram for clarity). We assume 
r^r = 10 -3 and show the histogram at four equally 
spaced time points, t/r = 5, 10, 15, 20. At short in¬ 
tegration times the Gaussian noise dominates, and the 
jump events contribute only a small bump distortion in 
the tail of the Gaussian. However, as the integration 
time increases, this bump grows to a flattened tail that 
becomes significant compared to the otherwise shrink¬ 
ing Gaussian noise. The overlap between P(7 | 10) and 
P(/ | 00) produces measurement error, which we discuss 
next. 


2. Measurement error probability 


We assume that the states |00) and 110) are discrim¬ 
inated by integrating the normalized quadrature I{t) 
and comparing the result I with a threshold value Ah- 
Slightly changing the notations used in the main text, 
we introduce the error probabilities for the two initial 


states as 


/»oo 

Perr(t|00)= / P(7 | 00) d7, 

7 7 t h 

(A.19) 

pith 


P err (f|l0)=/ P(7 10) d7, 

J —oo 

(A.20) 

where in the notations of the main text P err (^ | 00) 
and P err (t 1 10) = Per) , so that the overall error is 

— p(0) 

Perr (t) = ^ [P e rr{t \ 00) + P err (t | TO)] . 

(A.21) 


The error for the ground state 100) has a simple form, 


Perr (A 00) 


1 

2 


1 - erf 


1 + A h 

V 2r / t 


(A.22) 


which follows from Eq. (A.6). The error for the state 110) 
can be calculated as a weighted sum of contributions from 
zero, one, and two jumps, following Eq. (A.7), 

p err (t i io)= P0 p£\tm+ P1 p£(tm+P2P<g(t ito), 

(A.23) 

with the probabilities pk of k jumps given by Eqs. (A. 8) 
and (A.9), and the partial error probabilities, 

P f) (t | TO) = f ^ P (fc) (71 TO) dl, (A.24) 

.7—00 

obtained from the partial histograms P^ (7110) dis¬ 
cussed above. 

The zero-jump error is then similar to the ground-state 
error, 


Pen (t I 10) 


1 

2 


1 — erf 


cos(26>) - Z th 

V 2 T /t 


(A.25) 


and similarly to Eq. (A.22), it steadily decreases with in¬ 
creasing t. The exact formula for the single-jump error 
Pcrr (t\ 10) is very lengthy, but for practical purposes it 
can be significantly simplified. First, let us note that in 
the case T^. = F+ v , the distribution (A.14) for pW(7|10) 
is symmetric (since it is a convolution of a symmetric 
distribution Pi' j (z to t|10) and Gaussian noise). There¬ 
fore, for the symmetric threshold, Ah = 0, we have 
Porr (f 110) = 1/2. Second, for t/r > 4, the distribution 
P (1) (7|10) ~ [2cos(20)] _1 is practically flat for J t h close 
to zero, cos(2 9) — |/ t h| > \firji. Therefore, in this case 
dPirV(f|10)/d/ t h ~ [2cos(20)] _1 , which gives 

(A - 26) 

Very near the symmetric threshold, |/ t h| < r/3 1, a better 
approximation is possible, 


P e W(f|10) 


1 Ah cos(20) 

2 + 2cos(20) C1 v /277 V 


(A.27) 
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which corresponds to the exact derivative at It h = 0, 
P^(I = 0| 10) = [2 cos(20)]~ 1 erf[cos(20)/y / 2r/t]. Note 
that very good accuracy for pjrr (t\ 10) is not really 
needed since its weight pi in Eq. (A. 23) is small. In the 
case r s " w ^ T.+ , we can use approximation (A. 15) and ne¬ 
glect the Gaussian averaging, assuming cos(20) — |/ t h| > 
^/4 t jt. Then Eq. (A. 26) generalizes as 


P&\t I io) 


1 i th (r+ - r~ w )t ( _ il \ 

2 2cos(20) 8 \ cos 2 (20)/' 

(A.28) 


Note that for |r+, — -C 1 and I t h = 0 we have 
PeVr (t|10) ~ 1/2, which stems from the property that 
a nearly symmetric distribution (A.15) for remains 
nearly symmetric after convolution with the Gaussian 
noise. This explain the factor 1/2 in Eq. (38) of the 
main text, which follows from Eq. (A.23) with pq ~ 1, 
Pi ss T ~ w t, andp2 ~ 0. Similar approximations have been 
used in Eq. (43), in which we also assumed cos(20) « 1. 

The double-jump contribution to the total error (A.23) 
is very small because of small probability p 2 - Therefore, 


(O') _ 

it is sufficient to use a crude estimate for P e Vr (t|10). In 
particular, using Eq. (A. 17) and assuming p( 2 )(/|10) ~ 

Pi 2 ^(J|10), we find 


Pi r 2 r >(t|10) 


1 . Jth / t h 

4 2cos(20) 4 cos 2 3 4 5 6 (20) 


(A.29) 


Even though we neglect the three-jump processes, we 
actually take into account the main contribution from 
them automatically, by combining their probability p$ 
with pi in Eqs. (A. 8) and (A.9). This is because 
Pcrr(t|10) ~ 1/2 for I t h = 0 (following from the sym¬ 


metry of three-jump processes), which is the same as 
Porr (t|10). However, the dependence on I t h for the one- 
jump and three-jump terms is different. 

The blue dot-dashed line in Fig. 8 in the main text 
shows the error P err (i|10) calculated using Eq. (A.23) 
with the probabilities pk given by Eq. (A.8), the term 
Poi°r(t|10) given by Eq. (A.25), the term PirV(t|10) 
equal to 1/2 (because in Fig. 8 we use symmetric thresh¬ 
old, J t h = 0), and the term Per) (t | 10) obtained by inte¬ 
gration of the histogram (A. 18). This analytics fits the 
numerical result (red solid line in Fig. 8) significantly 
better than the simple analytics (dashed green line) dis¬ 
cussed in the main text. Actually, for the two-jump 
processes it is sufficient to use P e V/(t|10) = 1/4 [see 
Eq. (A.29)] instead of numerical integration; the result 
is almost indistinguishable. Note that this is practically 
equivalent to using p 0 = Pl = T~ w t [l-(3/4)T“ w f], 

P2 = 0. 

Finally, we emphasize that we have used the initial 
eigenstate 110) in the definition for measurement error 
in Eq. (A.20), since this is the optimal choice of logical 
encoding for the regime with T m -C A. If instead we 
use the bare state 110), then it will additionally collapse 
to a mixture of the eigenstates 110) and |01), which will 
increase the error, 

P err (t 110) = cos 2 (0) P err (f | ID) + sin 2 (0)[l - P err (f | TO)], 
= cos(20) P err (f | 10) + sin 2 6. (A.30) 

Thus, for the total error (A.21), we will have a nearly 
constant amount of additional error when distinguishing 
bare qubit states, P er r,bare ~ Pen + (p/A) 2 /2 for g « A. 
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